Particle-hole symmetric localization in two dimensions 
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We revisit two-dimensional particle-hole symmetric sublattice localization problem, focusing on the 
origin of the observed singularities in the density of states p(E) at the band center E = Q. The most 
general such system [R. Gade, Nucl. Phys. B 398, 499 (1993)] exhibits critical behavior and has 
p{E) that diverges stronger than any integrable power-law, while the special random vector potential 
model of Ludwig et al. [Phys. Rev. B 50, 7526 (1994)] has instead a power-law density of states with 
a continuously varying dynamical exponent. We show that the latter model undergoes a dynamical 
transition with increasing disorder — this transition is a counterpart of the static transition known to 
occur in this system; in the strong-disorder regime, we identify the low-energy states of this model 
with the local extrema of the defining two-dimensional Gaussian random surface. Furthermore, 
combining this "surface fluctuation" mechanism with a renormalization group treatment of a related 
vortex glass problem leads us to argue that the asymptotic low E behavior of the density of states 
in the general case is p(E) ~ E~ 1 e~' lnE ' , different from earlier prediction of Gade. We also 
study the localized phases of such particle-hole symmetric systems and identify a Griffiths "string" 
mechanism that generates singular power-law contributions to the low-energy density of states in 
this case. 



I. INTRODUCTION 

This article is devoted to a careful study of the localiza- 
tion properties of certain two dimensional systems with 
some special symmetry properties that make them stand 
apart from the other more generic universality classes 
of Anderson localization in two dimensions. The mod- 
els we consider are simply tight-binding models of non- 
interacting particles moving on bipartite lattices in two 
dimensions (2D). El The hopping matrix elements are ran- 
dom, but there are no on-site potentials, so the spec- 
trum is strictly particle-hole symmetric. Such 2DJ.C|cal- 
ization problems have received particular attentionBul in 
the context of the integer quantum Hall plateau transi- 
tion studies, and exhibit critical delocalized states only 
at the special energy E = 0. At all nonzero energies the 
eigenstates are exponentially localized with a finite local- 
ization length. 

Unlike the quantum Hall systems, which have an ana- 
lytic density of states (DOS) at the corresponding delo- 
calized critical points in the spectrum, the bipartite ran- 
dom hopping (BPRH) models have a strongly singular 
DOS at the band center. This unusual feature makes the 
sublattice problems very interesting in their own right 
since the disorder evidently has a very dramatic effect in 
the presence of the particle-hole symmetry. 

Recent studies of quasiparticle localization in dirty 
superconductors, where a somewhat analogous particle- 
hole symmetry appears naturally in the Bogoliubov-de 
Gennes formulation, have resulted in renewed interest 
in such special localization problems. For instance, the 
case with broken spin-rotation invariancea can be formu- 
lated as a pure imaginary random hopping (ImRH) prob- 
lem, and such ImRH problems constitute the simplest 
particle-hole symmetric generalization of the real BPRH 



models, obtained from the corresponding real BPRH by 
allowing additional purely imaginary hopping amplitudes 
between the sites of the same sublattice and performing 
a gauge transformation on one sublattice to make the 
inter-sublattice couplings imaginary. In the generic case, 
such systems all show some enhancement in the quasi- 
particle density of states at the band center, while the 
special case of a spinless superconductor with time re- 
versal invariance (which constrains the couplings within 
a sublattice to be zero) maps precisely onto a real bipar- 
tite localization problem with all the concomitant strong 
spectral singularities. 

To put our analysis of the bipartite models in context, 
it is useful to begin with a quick summary of what is 
known in the literature: Ludwig et ala showed that a 
special case of the sublatti ce pr oblem, the random vector 
potential model (see Sec. II 2 below), has a power-law 
density of low (i.e., near zero) energy states, 



p(E) ~ E- 



-l+2/z 



(1) 



parameterized by a continuously varying dynamical ex- 
ponent z. If we define an energy-dependent length scale 
L(E) from the density of states so that the integrated 
density N(E) = J^ dE' p(E') = L-' 2 (E), then Eq. (@) 
corresponds to dynamical scaling of the form E ~ L~ z . 
In related work, Gade and WegneiurEl developed a field 
theoretic description of the general sublattice problem 
(with random mass terms in addition to random vec- 
tor potential disorder — see Sec. |II 2| below), predicting a 
strongly divergent density of states 

p(E)-^e-^^\ (2) 

with x = 2, stronger than any finite- z power law Eq. (Fy) 
[these results have also been rederived very recently^ from 



a somewhat different field theoretic analysis] . This "Gade 
form" corresponds to a kind of activated (infinite-,?) dy- 
namical scaling 



\]nE\ ~ (InL) 



(3) 



when x > 1 . 

Finally, it is useful to note that these bipartite prob- 
lems are-jcipsely related to a well-studied vortex glass 
problemJZria 

The field-theoretic analyses that lead to these predic- 
tions, however, do not provide a direct physical picture 
of the origin of the divergent DOS at the band center and 
the nature of the low energy states. Developing such a 
picture is the main thrust of the present work. In partic- 
ular, our analysis leads us to suggest that the asymptotic 
scaling between the energy and length in the general sub- 
lattice problem has the infinite-.? form Eqs. (0) and (0), 
but with a different exponent 



X= 2 



(4) 



In the remainder of this section, we now outline our 
physical picture of the low-energy states, and summarize 
the basic argument that leads us to Eq. (||); for details 
and some relevant definitions, see Sees. Q and III below. 

The mechanism by which disorder produces a pile-up 
of low-energy states is best understood by first_consid- 
ering the special random vector potential modelo In the 
lattice version of this model, the hopping matrix elements 
are not independent random variables. Instead, they are 
given by the differences between the values of a random 
Gaussian field $(r) [see defining Eqs. (|14j) and (|15|)] on 
the two sites. This field has the statistics of a rough 
surface, with the variance of the difference in $ between 
adjacent sites being set by a disorder strength parameter 
g. As we argue later in Sec. Ill, the low energy states 
at strong randomness j> 1 "live" near the extrema of 
the surface $(r), and the logarithms of their energies 
are roughly given by the corresponding relative surface 
heights. Such extremal properties of two-dimensional 
Gaussian surfaces are well-characterized (see Appendix |A| 
for a summary of the relevant results); in particular, in 
a sample of size L, the corresponding prediction for the 
energy of the lowest state is 



I In El 



$max(£) - *min(£) ~ y/H In L , 



(5) 



implying power-law dynamical scaling with the dynami- 
cal exponent z growing as 



Vd 



(6) 



for strong randomness g^> 1. As mentioned earlier, Ref. [| 
indeed predicts a power-law DOS in this model, but with 
the dynamical exponent 



z = 1 H — 

7T 



(7) 



within the field theoretic analysis, this result appears 
to be perturbatively exact to all orders in g. However, 
the above "surface argument" shows that Eq. (p can- 
not hold for strong disorder, since the smallest energy in 
the problem cannot fall below the limits fixed by th e ex- 
trema of the surface. Indeed, we argue in Sec. IV A that 
there is actually a dynamical transition at g = g c = 2n; 
at this transition, the exponent z changes its behavior 
from the weak-disorder form (ffl) to the strong-disorder 
form (g) [see also Eq. (30)]. This dynamical transition 
is a counterpart of a static transitionHEil known to oc- 
cur when the wave function e~[*( r -'~* min ) (i.e., the zero- 
energy "pseudo-eigenstate" defined so that its peak value 
is 1) becomes normalizable in the limit L — > oo. 

Turning to the general case where the hopping ma- 
trix elements are independently random, we thus ex- 
pect some power-law density of states to be produced 
as long as there is some "vector potential" component 
to the randomness, as there is. Moreover, a field theo- 
retic renormalization group (RG) analysisQtl shows that 
any amount of a more general "random mass" disorder 
(see Sec. p 1| below) generates additional vector poten- 
tial randomness, driving the latter to stronger disorder 
at larger length scales and lower energy scales: 



3 cff (L) ~lni . 



(8) 



The nature of the low-energy spectrum and the corre- 
sponding dynamical scaling in the general case can now 
be heuristically understood by considering an "effective" 
random vector potential model with a scale dependent 
diso rder s trength parameter g c g(L) given by Eq. (0) [see 
Sec. IV B for a more precise argument from a different 
perspective using the results of Ref. 12 for a related vor- 
tex glass problem]. The low-energy states in this effec- 
tive problem are again associated with the extrema of 
some effective surface (which may be identified with the 
logarithm of the magnitude of the lowest-energy wave 
function), but the height extrema of this surface, and 
correspondingly the log-energies, now scale with L as 



| In SI ~ z e e(L)]nL ~ y/ gc g(L) lni ~ (InL) 3 / 2 , (9) 

which is the proposed "modified-Gade form" . [Note that 
the original Gade form would obtain if the weak-disorder 
form for z, Eq. (|7|), were used: 



!?! 



ln£| - z cS {L)\nL-rl- 5cff (L)lnL~ (lni)' 



(10) 



However, since g c s{L) grows indefinitely with L, this 
form cannot be used to obtain the true asymptotic be- 
havior.] 

So far, we have only alluded to critical bipartite sys- 
tems whose zero-energy states are not exponentially lo- 
calized. However, it is possible to drive such a system 
into a different localized phase keeping the particle-hole 
symmetry intact: This can be achieved, for example, 
by making some prescribed bonds that produce a com- 
plete dimer cover of the lattice stronger on average, thus 



introducing a preferred "dimerization" pattern. If this 
dimerization is weak, the system may remain critical, 
but strong dimerization will eventually drive it into a 
localized (insulating) state. We have found an inter- 
esting mechanism whereby disorder generates power-law 
contributions to the low-energy DOS in such bipartite 
insulating phases: The low-energy states are associated 
with the end-points of "strings" along which the back- 
ground dimer pattern is broken as in Fig. [H Estimating 
the relevant probabilities shows that this mechanism is 
actually operative for BPRH models in any dimension, 
as well as in the localized phases of the related ImRH 
problems (in the latter case, particular dimerization pat- 
terns appear quite naturally when describing disordered 
superconductors with broken spin-rotation invariance&3) . 
In our bipartite systems, this "Griffiths-like" mechanism 
can produce a power-law divergent DOS and dominate 
over all other mechanisms of filling the band gap near 
zero energy in such localized phases. 

Having introduced our principal results, we now con- 
clude this section with an outline of the rest of the pa- 
per: Sec. In] defines our models and reviews the connec- 



[II introduces 



tion with the vortex-glass problem. Sec. 
a strong-randomness picture of the low-energy physics, 
motivating the more precise bounds of Sec. IV on the dy- 
namical scaling in the system. Sec. |V| summarizes our nu- 
merical evidence in support of our analytical arguments. 
Finally, Sec. VI touches upon some unresolved questions 
and prospects for future study. 



II. MODELS, DEFINITIONS, AND THE 
CONNECTION WITH DIMERS 



A bipartite random hopping problem is completely 
specified by a single-particle Hamiltonian 



H = 



<a/3) 



t a p\a) 



H.c. = 







AB 



tAB 





(11) 



where a and (5 belong to sublattices A and B respectively 
of some bipartite lattice. In this article, we only consider 
real hopping problems in two dimensions, focusing atten- 
tion on the following specific models: ir-flux, honeycomb, 
and general bipartite lattice models. 

Each of these models is defined as a nearest-neighbor 
hopping problem on an appropriate lattice, with some 
additional constraints on—the allowed signs of the cou- 
plings: The 7r-flux modcL-3 is defined on a square lattice 
with a requirement that there is one half of a magnetic 
flux quantum per square plaquette. A convenient gauge 



jj+x 



(-l)^(j) 



'Tj+y 



ty(j) , (12) 



where j = {j x , j y } labels sites of the square lattice, x and 
y are unit lattice vectors, and t x (j) and t y { j) are nonneg- 
ative hopping amplitudes. The honeycomb (brickwall) 



lattice modelE-3 is defined on a honeycomb lattice with 
nonnegative hopping amplitudes. Finally, by the general 
lattice model, we mean a nearest-neighbor hopping model 
on a square or honeycomb lattice with no constraint on 
the allowed signs. 

Naturally, the pure-system spectra in the above models 
can be quite different: The pure 7r-flux and honeycomb 
lattice models have Dirac points and a linearly vanishing 
density of states at E — 0, while the pure rectangular lat- 
tice model has a Fermi surface with a smooth DOS near 
E = in the general anisotropic case and a van Hove 
singularity in the isotropic square lattice case. However, 
all these models are believed to have similar low-energy 
long-wavelength localization physics in the presence of 
randomness. In our studies, it is often convenient to 
work with a particular model when using a particular ap- 
proach, and we switch among the models a fair amount 
in what follows. 



1. Continuum description of models 

The pure 7r-flux model has two Dirac points. The 
corresponding weakly-disordered model is described by 
a continuum Dirac Hamiltonian, which may be conve- 

hi 



niently written asEJ hi 



hi 



with 



h2 = a x (-iv x d x +iA x ) + o~ y (-iv y d y + iA y ) - iV + Ma z . 
In the above, the slowly varying vector and scalar po- 



tentials are implicitly defined by 2t 



(0) 



2t 



(<>) 



v v , -2St x (j) = MjK-i) 3 " +3v + V(j)(-1)" + .., 

-26t y (i) = A y (i)(-iy*+Jv +M(i)(-iyy +..., and the 
lattice constant has been set to unity. A similar con- 
tinuum Hamiltonian obtains for the honeycomb lattice 
model. 

The Dirac Hamiltonian h 2 is non-hermitian and con- 
tains a real random mass M, a random imaginary po- 
tential — iV , and a random imaginary gauge field. The 
{M, V} terms form a "complex random mass" part, 
which we refer to simply as the random mass part of 
the disorder. In its absence, M = V = 0, the full Her- 
mitian Hamiltonian /14 decomposes into two blocks, with 
the "1-4" block having a form 



h-A = o- x {-id x + Ai) + a y {-id y + A 2 ) 



(13) 



identical to the (real) random vector potential model of 
Rcf. H [here {Ai,A2} = {A y , — A x }], and similarly for the 
"2-3" block. Here and henceforth, we will therefore refer 
to the {A Xl A y } part of /14 as the random vector potential 
part of the disorder. 



2. The random vector potential model on a lattice 

Consider a "random-surface" BPRH model con- 
structed from a field $(r) and a pure system Hamiltonian 
with bare hopping amplitudes t\^ as follows: 



ta/3 



P *(o) f (°) e -*(/3). 



v a(3 



(14) 



J(o) 



One can easily see that if ^ A ' is a zero-energy 
eigenstate of the pure model, then ^ a with ^ a(o) — 
e~*( a 'Hf\' (a) is a zero-energy eigenstate of the random- 
surface model. Moreover, it is easy to see that this con- 
struction provides a particular lattice realization of the 
random vector potential model, if the pure system one 
starts with is the non-random 7r-flux model on a square 
lattice. Indeed, in the continuum limit, we obtain M = 
V = 0, A x = v x d x $>(r), and A y =v y dy$>(r), so that h^ of 
Eq. jl^ ) becomes h_A — o- x (—id x +d y $>)+o-y(—id y —d x $). 
Our random surface $(r) thus represents the physical 
(non-gauge) degrees of freedom of the random vector po- 
tential, and we choose it to be Gaussian 



Prob[$] oc exp 



■% !**&*? 



(15) 



with dimcnsionless disorder strength g\ the relevant prop- 
erties of such surfaces in two dimensions are summarized 
in Appendix [Al 

A remark is in order here: In the continuum limit of the 
general problem, we have seen that there is a clean sep- 
aration of the disorder into its random vector potential 
and random mass parts. However, in the original model 
on a lattice, such a separation is not at all obvious — the 
only precise statement we make at the lattice level is that 
a model has purely vector potential randomness when it 
can be put in the form given by Eq. (|14|). Renormal- 
izing on a lattice in Sec. Ill, we will still speak of the 
flows of the two parts of the disorder, implicitly assum- 
ing some coarse-grained description and relying on the 
corresponding results for the continuum model. A more 
precise lattice alternative to such a hybrid RG approach 
is the domain wall energetics picture of Ref. [La which 



we will have occasion to appeal to in Sec. IV B 



where the sum is over all dimers in this covering, and 
e Q ^'s are the random bond energies. The fermion-dimcr 
connection is stated most easily for open boundary con- 
ditions: In this case, the partition function of the dimcr 
model Zd = J2m ex P ( — £<2pVt]/Ty at a given dimer tem- 
perature Td can be written as a determinant of an appro- 
priate connection matrix tAB, Zd[A, B] = det tAB, where 



t a p = ±e" 



„/3 /T d 



(17) 



and the signs of the hopping amalitudes are chosen ex- 
actly as in the 7r-fiux model.Ej'tH The honeycomb lat- 
tice model with nonnegative bonds has a similar dimer 
connection (our definitions of the 7r-flux and honeycomb 
models were made precisely with this in mind), but the 
more general lattice model has no such direct connection. 



4- Localized bipartite hopping systems 

The commonly studied bipartite hopping models (e.g., 
all of the above models with some generic distributions 
of random hopping amplitudes) are found to be critical, 
by which we mean that the behavior at E = is nei- 
ther truly delocalized nor truly (exponentially) localized. 
Abusing language somewhat, we will often refer to a sys- 
tem in such a critical state as being in a "delocalized" or 
"metallic" phase — this serves to emphasize the distinc- 
tion between such a phase and truly localized insulating 
phases that obtain, e.g., by introducing strong dimeriza- 
tion in the couplings. The pure-system counterparts of 
these localized phases are typically gapped band insula- 
tors. In our studies, we will specifically consider only 
one such localized phase produced by introducing "stag- 
gered" dimerization on the square lattice (see Fig. |l]) ; the 
honeycomb lattice version of this is shown in Fig. 13(a). 



III. ON THE ORIGIN OF LOW-ENERGY STATES 

We begin by noting that the Gade-like dynamical 
scaling Eq. (g) with x > 1 suggests that the effec- 
tive value of disorder scales to infinity at asymptotically 
low energies (dynamical exponent z = oo), albeit very 
weakly. This nieiisjates us to try a strong-randomness 
RG descriptionMlB 



3. Connection with random dimer problems 



There is a direct connection between the 7r-flux random 
hopping problem and a random dimer model on the same 
square lattice. The dimer model consists of all complete 
coverings (matchings) {M} of the (bipartite) lattice by 
the nearest- neighbor dimers; the energy of a given cover- 
ing JVC is 



e d [M]= J2 e ^< 

{ap)eM 



(16) 



A. Strong-randomness RG 

Consider the bipartite hopping Hamiltonian (|ll|). The 
eigenstates of H occur in pairs with energies ±25, and 
the strong-randomness RG proceeds by eliminating, at 
each step, such a pair of states with energies near the 
top and the bottom of the band: One finds the largest 
(in absolute value) coupling in the system, say £1.2 con- 
necting sites IS A and 2g_B; this defines the bandwidth 



2fi = 2max{|i Qj ^|}. If the distribution of the couplings is 
broad, the two eigenstates of the two-site problem H[l, 2] 
with energies ±il will be good approximations to eigen- 
states of the full H, since the couplings t\^ and i2 jQ . of 
the pair to the rest of the system will typically be much 
smaller. These couplings can then be treated perturba- 
tively, and eliminating the two high-energy states living 
on the sites 1 and 2 gives us the following effective cou- 
plings between the remaining sites: 



^q : /3 — ta,/3 



t n 



*i 2 *1, 







(18) 



The renormalized Hamiltonian is again a bipartite hop- 
ping problem, but with two fewer sites; in particular, no 
diagonal terms are generated. This procedure is then 
iterated to reach lower and lower energy. 

As a first attempt, we implemented this RG numer- 
ically in both localized and metallic phases. In the lo- 
calized phases, the RG indeed provides a qualitatively 
and quantitatively accurate description of the low-energy 
physics, with the relative "decimation errors" diminish- 
ing quickly at low energies; in this case, the RG is thus 
a consistent scheme that correctly describes the Grif- 
fiths effects analyzed in detail below. However, in the 
delocalized phase, the observed flow to stronger disor- 
der is very weak, and the consistency of the RG is less 
clear. We do not address in detail here this question of 
self-consistency (that is, in which regimes — if at all — the 
system flows to stronger disorder, with the RG becom- 
ing more and more accurate). Rather, in the metallic 
phase, we view the RG as providing a heuristic picture 
that complements the more precise arguments we present 
later. Note, however, that the RG actually contains some 
exact .and useful information for any strength of the dis- 
order:^ The rule Eq. ( J18| ) is an exact transformation for 
the zero-energy wave function, and as such provides use- 
ful information on the E=0 localization properties. Also, 
running this RG corresponds to employing the Sturm se- 
quence method for calculating the integrated density of 
states at E = using a particular order that is, with- 
out any a priori knowledge, numerically most reliable; 
the intermediate terms that appear in this process at a 
particular length scale give us a rough idea about the 
corresponding energy scale in the problem. 



B. On the origin of states 

The strong-randomness RG associates the low-energy 
states of the original Hamiltonian with the sites that 
"survive" down to the corresponding energy. Such a 
"free" site — say, an A-s\te — is found if all the -B-sites in 
its immediate neighborhood happen to be "locked" (by 
stronger bonds) into pairs with some other ^4-sites. Then, 
this free A-site can only couple to the next available B- 
sites, which are far away, and the corresponding effective 
couplings are relatively weak since they are mediated by 
a substantial intervening pair-locked region. 
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FIG. 1. "String" mechanism for generating low-energy 
states in the band insulator phases. Filled and open circles 
represent A and B sites, respectively. Medium thick bonds 
represent a particular background dimer pattern (i.e., bonds 
that are typically stronger), broken bonds from this set rep- 
resent couplings that happened to be atypically weak in a 
given disorder realization, while the heavy thick bonds repre- 
sent atypically strong couplings. In the above figure, there is 
a pair of low-energy states associated with the end-points of 
the resulting string. 



It is at first difficult to see how such free sites and the 
corresponding low-energy states can be produced often 
enough in dimensions d > 1: Naively, the energy scale 
associated with a pair {a, /?} of "isolated" sites distance 
I apart is ~e~ CEl . Now, it may seem that to produce such 
a pair we need of order l d specific events — for example, 
we can imagine having a particular pairing pattern for all 
the other sites in the neighborhood of the two unpaired 
sites a and (3 — and this has a very low probability Prob ~ 






in d > 1 . 



In the dimerized band insulator phase, however, there 
is an insulating background to start with, and to produce 
such a pair one only needs to "break" the background 
dimer pattern along a string joining the two sites — i.e., 
of order I specific events (see Fig. nl), with the resulting 
occurrence rate ^e c ? high enough to give a power-law 
contribution to the low-energy density of states. Our 
numerical RG and exact diagonalization studies confirm 
that this is indeed what happens in such band insulator 
phases. Furthermore, this "Griffiths-string" argument 
goes through for localized phases in any dimension. 

In the "delocalized" phase, on the other hand, we do 
not have such an insulating background. However, the lo- 
cal properties of the "isolating" region need not be very 
"sharp" , as was assumed in making the estimate for the 
probability of occurrence Prob ~ e~ Cpl . Rather, there 
might be some smooth "textures" of the random pair- 
ing ( "dimer" ) pattern producing such isolated sites ( "de- 
fects"), and these may occur with much higher probabil- 
ity. 

To be more specific, first consider the "random- 
surface" model Eq. (HJ) (i.e., the random vector poten- 
tial model) in the limit of strong disorder g> 1. The 
RG, being exact for the zero-energy wave function, pre- 



serves the random-surface character of the problem leav- 
ing $(r) unchanged: the randomness in the effective cou- 
plings that are generated by the RG rule Eq. jlq ) is of the 
same random-surface form Eq. ( fA\ j with the same $(r). 
Thus, at each stage of the RG, the magnitude of the cou- 
pling between any two remaining sites a and (3 is roughly 
^ e *( Q )-*(^)_ Since the RG eliminates the stronger cou- 
plings first, an A-site a m i n that is a local minimum of 
<f>(r) cannot be decimated out until the "probing" log- 
energy scale is large enough to "see" the next even lower 
minima: Indeed, before this happens, the other >l-sites 
in the neighborhood of a m i n , having higher $, are more 
strongly coupled to any available -B-sites in the same 
neighborhood. The same holds for any i3-site /3 max that 
is a local maximum of $(r). As a result, within our 
strong-randomness RG picture, the low-energy states of 
the problem are associated with the local extrema of the 
surface $(r), while their effective log-energies are given 
by the corresponding surface height differences. In par- 
ticular, in a finite system, the smallest energy scale is 
associated with the global extrema of <I>(r). The argu- 
ments leading to Eqs. (|])-(@) °f the Introduction now 
follow. 

Thus, within our somewhat heuristic RG approach, 
"isolated" sites are produced often enough to result in 
a power-law contribution to the density of states, with 
the dynamical exponent z ~ ^/g in the strong random- 
ness regime. Here, we are ignoring effect of running the 
RG transformation on the pure model itself — this will 
change the actual value of the exponent z, but not the 
~ y/g disorder contribution for 9 > 1. For more pre- 
cise bounds that complement these heuristics and lead 
to identical conclusions, see Sec. IVA. 



In the general case, the heuristic RG arguments 
Eqs. (pi)— (@) in the Introduction now follow and yield 
the modified-Gade forms Eq. (g) and (||) with x = 3/2; 
more precise argum ents l eading to identical conclusions 
are presented in Sec IV B. However, it should be empha- 
sized again that we have not found a clear-cut separation 
of the random vector potential and random mass parts 
of the disorder at the lattice level. We therefore have 
to rely £>n the results of a different field-theoretic RG 
analysisLrti when discussing how the two types of ran- 
domness flow and affect each other (note that a similar 
somewhat heuristic argument was used in Ref . [19 to ob- 
tain crossover scales in the vortex glass problem). 

We conclude this discussion of the origin of low-energy 
states with some miscellaneous observations: Recall that 
the surface in the random vector potential case is also 
the logarithm of the pseudo-zero-energy-wave function 



*aO) 



-*( Q ) . Given that the RG rule is an exact 



transformation for such a wave function in the general 
case as well, the above analysis suggests that the E 7^ 
low-energy states in the general problem are again re- 
lated to the extremal properties of this E = wave func- 
tion. It is also useful to note the predictions of this "sur- 
face" mechanism in other dimensions. In one dimension, 
any disorder is of purely "longitudinal" (random-surface) 



type, and the Gaussian surface is more rough, with 

I ln£| ~ $ max (i) - $ min (L) ~ VI ; (19) 

this is precisely the Dyson critical scaling for the ID 
chain. For bipartite systems .that are finite in all but one 
direction (such as ladderaH3 : EJ) , the DOS again has the 
Dyson form at special delocalized critical points in the 
phase diagram, since the corresponding surface "looks" , 
at long wavelengths, like a ID Gaussian surface. On the 
other hand, in dimensions three and higher the Gaus- 
sian surface is basically flat, and this mechanism can- 
not generate power-law contributions to the density of 
states. Much of what has been said above for the two- 
dimensional problem is therefore very specific to 2D. 



C. From low-energy states to optimized defects 

We now establish a suggestive connection between the 
low-energy states in the bipartite hopping problem and 
monomers (defects) in the corresponding random dimer 
problem. The dimer connection leads us to a heuristic 
"dimer RG" prescription for the low-energy states — this 
dimer RG approach provides us with a very useful phys- 
ical picture that underlies the more precise analysis of 
Section £y|. 

Consider running the above "zero-energy" RG in some 
specified order. (The naive strong-randomness RG de- 
scribed earlier uniquely prescribes a particular order by 
picking the strongest effective bond at each step, but it 
is important to realize that in the absence of any rapid 
flows to stronger effective disorder, there is nothing "sa- 
cred" about this ordering, and other sensible choices are 
possible.) We denote the A- and -B-sites decimated out 
at a given stage of the RG as 21 = {01, . . . , a n } C A and 
05 = {61, . . . , b n } C B respectively. The effective cou- 
pling between any two remaining sites a and (3 is given 
by 



t a , p — ta,fl — t a , <8 *<a, B ^21. P 



(20) 



Here, t^^ is a matrix of the original bonds joining two 
subsets <t C A and £> C B (thus, t a ^ is a row-vector, 
t%,p is a column- vector, and £21, s is an n x n matrix). 
Expression (J2(j) follows readily if we note how the RG 
"solves" for the zero-energy wave function by eliminating 
any reference to the decimated sites, and this holds for 
any order of the decimations. Furthermore, it is easy to 
check that this may be rewritten as 



t' 



det i a+ % y p + <B 



a, P 



det ia 



(21) 



<8 



where a+2t = {a, ai, . . . , a„}, and similarly for /3+Q3. It 
is this form that we find useful below. 



and S of equal cardinality 



Now, recall from Sec. II 3 that for any two such sets € 



detic.s = ^sgn(M) 

M 



11 ta P ' 

(o!/9)e3vt 



(22) sponding estimate 



where the sum is over all complete dimer coverings of 
£ U 55, and sgn(M) denotes the appropriate permuta- 
tion sign for matching M. Each term of the sum has a 
form ± exp ( — fi^pVCJ/T^J, with the dimer bond energies 

e a /3 defined from \t a p\ — f2oe~ eQ,3 ' Td . In the limit of infi- 
nite randomness, or equivalently, zero temperature of the 
dimer system Td = 0, one term — the ground state dimer 
covering — dominates the whole sum; obviously, the per- 
mutation signs play no role in this limit. 

Consider first this infinite-randomness limit. The 
largest term in det t a+ <& : 0+9$ corresponds to the optimal 
(ground state) dimer covering of {a+2l} U {/3+25}, while 
the largest term in det ia, *B corresponds to the optimal 
covering of 21U23. The two optimal coverings "differ" by a 
"string" connecting the "added" sites a and (3, while the 
effective coupling t' a « in this infinite-randomness limit, 
Eq. (|l]) , is determined precisely by "propagating" along 
the string t' a a = tit 3 . . . t2fe+i/(^2 ■ ■ • hk)- Note that this 
string is the optimal (i.e., corresponding to the largest 
possible \t'\) such path from a to (3 utilizing the sites 
decimated out. 

This suggests the following infinite-randomness "dimer 
RG" sequence: At the n-th stage of such a dimer RG 
analysis, we find the ground state configuration M„ 
with n dimers, defining £n = £<j[M„ ]. The 2n sites 
covered by dimers at this stage are to be thought of as 
having already been decimated. At each step of this 
RG, precisely two sites are "integrated out" since the 
optimal (n + l)-dimer covering differs from the optimal 
?i-dimer covering by two added sites (and the connect- 
ing string). We now need to provide, for the random 
hopping problem, some notion of the bandwidth (i.e the 
value of the cutoff energy) corresponding to each stage. 
A natural choice is to associate the effective hopping 
amplitude \t'\ = ft exp [ - (Sjf+1 - t^ s) )/T d ] with the 
two sites "eliminated" in going from stage n to n + 1 
this then specifies the value of the cutoff in the hopping 
problem at this RG step. The monotonicity property 
£„+i — £n > £n — £„_i guarantees that the cutoff in- 
deed decreases monotonically as the RG proceeds. Note, 
however, that the "pairings" of the sites decimated out 
are somewhat subtle; these do not stay rigid as the RG 
proceeds, and can change each time we decimate more 
sites. Furthermore, there is no simple relationship in 
general between the sequence of decimations within this 
dimer RG, and that in the naive strong-disorder RG se- 
quence described earlier. 

Of course, this can all be reformulated starting from 
the ground state complete dimer covering (with N = 
L 2 /2 dimers) of the whole lattice, and then adding pairs 
of monomers (defect pairs) into the system in the most 
optimal way. In particular, within this heuristic picture, 
the lowest energy state in the hopping problem corre- 
sponds precisely to the first such pair, with the corre- 



r 



n(dimcrRG) _ 



z~v / vopt.dcf . \ 



(23) 



here v pt.def. — &jy—i c*jv 



?(9 S ) 



Em"' is the optimized (lowest en- 



ergy) defect pair energy, and all of the foregoing assumes 
the system has complete coverings to start with£3 

Returning to the finite-randomness case, there are sev- 
eral issues that stand out when trying to relate the low- 
energy states of the hopping problem to the monomers 
in the corresponding finite-temperature random dimer 
system. The determinant det t%. <g and the dimer par- 
tition function Zd[%l, 23] still contain the same (in abso- 
lute value) terms, but now the dimer ground state by 
itself does not determine the corresponding sums, so the 
permutation signs become important. Even if we use 
the dimer partition function only to bound the absolute 
value of the determinant, we are still faced with the dimer 
problem at finite temperature, which now has significant 
entropic contributions, since local moves from the ground 
state cost only finite energy. In this situation, we can no 
longer think solely in terms of the ground state pairing 
when attempting to estimate ^[21, 23]. 

However, we are primarily interested in the ratios of de- 
terminants like det t^s and dettA- a , B-p, i-e., the dif- 
ference in the corresponding free energies of the dimer 
system without defects and with two defects. Since 
a large part of the entropic terms in the two systems 
comes from the same local moves, the ground state en- 
ergy difference of the two may still be a relevant ob- 
ject. Put another way, we need to consider the effec- 
tive thermodynamics of the defects themselves — i.e., with 
all free energies measured relative to the defect-free sys- 
tem. As we sh all see within a more precise formula- 
tion in Sec. IV B| , this effective thermodynamics is indeed 
ground-state dominated, allowing us to make fairly pre- 
cise statements for finite values of bare disorder as well. 



IV. "VARIATIONAL" BOUNDS ON THE STATES 
NEAR THE BAND CENTER 

Motivated by these heuristic RG considerations, we 
now establish more rigorous bounds on the dynamical 
scaling in the system. Consider the bipartite Hamilto- 
nian (pTj) in a finite sample and assume that there are no 
zero-energy states and no edge states; a concrete realiza- 
tion would be some regular lattice block with some nice 
boundary conditions, e.g., the 7r-flux model on an L x L 
square with free be and even L. The smallest positive 



energy E B 



E^ lin of H is then given as 



E n 



1 



11*^11 



(24) 



Here, ||C|| is the matrix norm of an operator C, and 
satisfies the inequalities 



max|C y |<||G||< [ZlCij? <J2\C i: 



(25) 



'■./ 



hJ 



Expressions (E4) and (pq) formalize our earlier, more 
heuristic association of low-energy states with optimized 
defects in the dimer problem. Indeed, we have 



(^4 B )/3a - ± 



dettA-a,B- 

det i A B 



(26) 



which is precisely the previously encountered ratio of the 
determinants, and the above inequalities suggest that the 
two "defect" sites a and /3 should be chosen in some 
optimal way. 

Below, we use these inequalities to obtain, most impor- 
tantly, lower bounds for the smallest positive energy.E!! 



A. The random vector potential model 

Applying the above bounds to the "random-surface" 
model Eq. (|14|), we obtain, in particular, 

||f^ 2 <£e 2 ^- 2 ^>|G (/3,a)| 2 , (27) 

a.,6 



where G (/3, a) = (i 1 ) / g Q 



|-ff _1 |a) is the E 





Green's function (propagator) for the pure lattice prob- 
lem. The propagator Go can be calculated explicitly for 
each particular pure system. In the following, we con- 
sider specifically the 7r-flux lattice model with free be; the 
lattice propagator (Dirac fermions) can be found, e.g., 
in Ref. |4], but here we need only note that Go(i"i,r 2 ) 
is bounded and behaves as |Go(ri,r 2 )| ~ |i"i — r 2 | _1 
for large |ri — r 2 | 3> 1. A simple (completely rigor- 
ous) bound Pxsll — cons t x 5(L) immediately follows 
from the boundedness of Go- Here, 3(L) is the parti- 
tion function for a classical particle living in a reduced 
potential V(v)/(kBT) = 2$(r) — see Appendix [K], where 
we also summarize other relevant properties of such 2D 
Gaussian surfaces. Furthermore, it is actually possible 
to provide stronger bounds — below we discuss this sep- 
arately for the strong-randomness (g > g c = 2ir) and 
weak-randomness (g < g c ) regimes (see Appendix W for 
the significance of g c and the distinction between the two 
regimes in the random-surface problem). 

In the strong-randomness regime, we expect that the 
sum (|27| ) is dominated by the /?'s near the maxima of 
$(r) and the a's near the minima of <&(r) [analogous to 
Eq. ( |A3| ) for 3(i)] ■ We also expect that these extrema are 
separated by distances of order L (the precise statement 
would be that the distance between the global maximum 
and minimum of ^(r) is L times a random number of 
order one). This suggests that 



^abII — (num.fact.) x 



,(i)-*min(£) 



3(i) 



(28) 



The expected uncertainty in the rhs of this expression 
from one disorder realization to another is simply a ran- 
dom O(l) numerical factor (num. fact. ). [Note that even 
though the above estimate of the sum (|27|) does not con- 
stitute a rigorous proof, it is fairly robust — moreover, a 
more formal proof is likely possible along the lines of 
Refs. |lO],[ll].] Since the sum in the estimate is essentially 
dominated by few terms, this also provides a lower bound 
for ||t^g||, that differs from the upper bound (Eq) only 
by an O(l) numerical factor: 



\\i A B\\ >max|(£ A s)/3c 



: const x 



K (i)-*mla(i) 



the bound is rigorous (and completely general) since the 
distance between the global extrema cannot exceed the 
sample size. Thus, the smallest positive energy is ex- 
pected to scale with the system size as 



E n 



(num. fact.) x 



L 



(lnL): 



3(i) WiF- 1 

In particular, for the dynamical exponent we obtain 



(29) 



9 > 9c 



(30) 



This is our central result; pictorially, the surface ^(r) 
simply does not allow the low-energy states to fall be- 
yond its extrema. [In Eq. ( |29| ) we have also included 
the logarithmic correction from Ref. O to show that the 
minimal energy is somewhat larger than expected from 
the simple-z power-law; note, however, that the actual 
form of such log-corrections caa, in general, depend on 
the boundary conditions used.lill] 

In the weak-randomness regime, the bound Eq. (|27j) is 
not as precise. For instance, in the pure-system {g = 0) 
case, we obtain 



J2\Go(P,a)f 



= TriT n - 



L 2 \uL 



(31) 



a,p 



i.e., i?min > L _1 (lnL) -1 / 2 , which gives the correct expo- 
nent z = 1, but is certainly an underestimate. Similarly, 
in the g =/= case, we can prove that for any r\ > 



U1b\ 



< const X Ll+2ff/9c+r, 



(32) 



with probability one in the limit L — ► oo. The proof 
(following the lines of Ref. 10) applies the inequality 
Prob(X > K) < (X)/K, valid for any nonnegative ran- 
dom variable X, to the particular X equal to the rhs 
of Eq. (J27]) . Thus, the dynamical exponent z obeys the 
following bound quite generally: 



z < 1 + — 
9c 



(33) 



Note that this result is valid in the strong disorder regime 
as well — however, in this case an individual realization es- 
sentially never samples the tails of the distributions that 




FIG. 2. Density of states in the 7r-fiux model with the "ran- 
dom-surface" disorder [see Eq. (UM] calculated for 192 x 192 
system with open boundary conditions, for varying disorder 
strength g. The data is averaged over 10 disorder realiza- 
tions. Note that the plots have some curvature for strong 
disorder — this is probably because of the logarithmic correc- 
tions to the power-laws in this regime. 



■ z (numerical DOS) 

- z(g<g c )=1+2g/g c 
z(g>g c )=4(g/g c ) 1/2 -1 

- z(g>g c ) with log-correction 
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FIG. 3. Dynamical exponent for the 7r-flux model with 
the random-surface disorder extracted from the data of 
Fig. @. The additional lines show theoretical predictions in 
the two regimes (with g c = 2w)), and also a rough esti- 
mate z c ff(L) — z(oo) — (3/2) y/g/g c (In In L)/ (In L) of the 
log-correction in the strong-disorder regime for the specific 
system size studied. 




determined the average (A), and our earlier treatment 
leading to Eq. ( p0| ) is more appropriate and provides a 
much sharper bound. 

For weak randomness g < g c , we can also give a plausi- 
ble argument for the opposite inequality z > l + 2g/g c . It 
suffices to provide an upper bound for -Emm, and a good 
trial wave function for this seems to be 



ipB = t AB e 



-*(A) 



with E m[n < WtABipBW/WipBW bounded by 
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(34) 



(35) 



To translate this into an actual lower bound on z, we first 
estimate the sum over a sites in the denominator as 



E G °(^ 



ft e 



-2*(a) 



E 



G (/3,a)| e - 2 *(«) . (36) 



This estimate is clearly reasonable, since Go( r i> r 2) is 
some fixed function, while the sums involving e -2 **-"' 
over macro-subregions (with Go > and Go < 0) always 
have O(l) random numerical factors [see, e.g., expres- 
sions for 3(L) in Appendix]. Now, using |Go| > L^ 1 we 
obtain E m { n < (num. fact.) x L/3(L), which gives the de- 
sired coinciding lower bound for z in the weak-disorder 
regime (note that this also provides a similar coinciding 
lower bound in the strong-disorder regime, and thus re- 
produces our result for z in this case). 

We conclude with two comments on the results ob- 
tained here. Firstly, although we cannot give a proof, 
we expect that the expression for z in the weak-disorder 
regime, 



25 
9c 



9 <9c 



(37) 



is exact without any logarithmic corrections. Our ex- 
pectation is bolstered by the fact that Eq. ( p7| ) coin- 
cides with the original result of the weak-disorder replica 
analysis of Ref. [which is known to give exact re- 
sults also for the static quantities of Appendix |X] in the 
regime g < g c ]. Moreover, the two expressions for z(g), 
Eqs. ( p0| ) and (f37|), join continuously and smoothly at g c 
[and we expect that E m i n — (num. fact. )x L/3(£) is valid 
for any disorder strength]. 

Secondly, note that the above results are, strictly 
speaking, valid only for the exponent z\ that determines 
the scaling of the smallest gap in the system. However, 
this is not expected to be a serious issue: Thus, having 
the upper bounds on i? m j n (L), we are guaranteed a bulk 
density of states with the DOS dynamical exponent zdos 
at least as large. Moreover, we expect quite generally 
that there is only one dynamical exponent in the system 
zrjos — z \ = z - Some numerical checks for the bulk den- 
sity of states are shown in Figs. || and [| and these are in 
good agreement with the above expectations. (See also 
Ref. |25| studying essentially identical dynamical freezing 
phenomenon for random walks in random log-correlated 
potentials in one and two dimensions.) 



B. Scaling of the lowest energy in the general case 

In general, we can always bound the numerator 
in Eq. (|26| ) by the corresponding dimer partition function 
with two fixed monomers a and (3 as | dettA- a .B-(}\ < 
Zd[A — a, B — (3\. We now specialize to the 7r-flux model 



with free be (or any other model with a direct connec- 
tion to dimers), so that for the lattice with no defects 
we have | deti AB \ = Z d [A, B] = Z N , with TV = L 2 /2 the 
total number of dimers. In this case, we do not have to 
worry about the denominator at all, obtaining, e.g., 



\t 



AB\ 



< 



J2 a0 Z d [A-a,B-f3] Zjv-i 



Z d [A,B] 



-N 



(38) 



where Z/v-i is the partition function of the dimcr prob- 
lem with two free defects allowed to be anywhere on the 
lattice. 

Now, this ratio of partition functions describes the ef- 
fective finite temperature thermodynamics of a pair of 
defects (i.e., with all contributions to the free energy 
measured relative to that of the dimer system with no 
defects) . The analysis of Ref. |lj suggests that this effec- 
tive finite-temperature thermodynamics is equivalent to 
the thermodynamics of a pair of classical charged parti- 
cles (with opposite charges for A and B monomers) in a 
random "electrostatic" potential v(r) with average cor- 
relations 



[v(r) - v{r')] 2 



ln 2 |r 



(39) 



The extremal properties of this random surface v(r) have 
also been characterized in Ref. O. In a finite box of linear 
size L, we have for the average and the standard deviation 
of the distribution of the global extrema 



tW ~ (lni)^, a(tw) ~ (lni) 1 ^ , (40) 

Reasoning as in Ref. [0], we see that these extremal prop- 
erties immediately imply 



-'JV-l 

Zjsr 



< L exp 



T d 



< exp 



const x (lnL) 3/2 
(41) 



and correspondingly for the lower bound on the smallest 
positive energy 



ln.E min > -const x (InL) 3 / 2 



(42) 



In the above, making the positive number "const" slightly 
larger takes care of any prefactors L r and any uncertainty 
o((lnL) 3 ' 2 ) in f max . This result is a simple manifestation 
of ground state dominance^ for the finite-temperature 
partition sum of particles moving in the random poten- 
tial v(r) with correlations Eq. (p9[). The above crude 
estimates are good enough to reach the main conclu- 
sion Eq. (|4|) precisely because of this, even though the 
ground state dominance is a very weak one. Since a 
few defect positions effectively determine the above es- 
timates, we expect that Eq. ( ff2| ) provides not only the 
lower bound but the actual leading term in the logarithm 
of the smallest positive energy in the system. This im- 
mediately leads to the modified-Gade scaling proposed in 
the Introduction — our analysis here thus provides strong 



evidence in support of the RG arguments outlined in the 
Introduction. 

To get some idea about the subleading terms and, in 
particular, the relationship between the actual lowest 
gap, the intermediate bound Eq. (38) obtained from the 
effective thermodynamics of a free defect pair, and the 
heuristic dimer RG "estimate" Eq. (p3J) from the opti- 
mized defect pair, it is instructive to consider, from this 
perspective, the random vector potential model in the 
strong-randomness regime g 3> 1 . 

Recall first our earlier direct estimate of ||tTg||) i.e., 
Using the Dirac fcrmion propagator 



the lhs of Eq. 
for £q , we have 



I(*abWI 



const 

l r /3 - r Q 



x e 



*(/3)-*(a) 



(43) 



Invoking the ground state dominance in the strong- 
randomness regime, we obtain ||i^ B || ~ x _1 e* max_ * min . 
Consider now the relative partition sum Zn-i/Zn for 
a free defect pair, i.e., the rhs of Eq. (p8|). In this 
"random-surface" case, for a given fixed configuration of 
monomers, any dimer covering has the same energy, and 
the partition function is essentially purely configurational 
(entropic). Thus, for two fixed defects a and /3, we have 



Z d [A-a,B-0\ 
Z d [A,B] 



const 



hs-r a |V2 



*G9)-*(a) 



(44) 



since the relative number of dimer configurations decays 
as |r^ — r Q | -1 / 2 (Ref. ^J) with increasing separation be- 
tween the monomers. Taking into account the ground 
state dominance in the strong-randomness regime, we 
then estimate Zn-i/Zn ~ £~ 1 ' 2 e* ma - x ~* min . Thus, the 
rhs of the inequality ( p8| ) overestimates the lhs by a factor 
of L 1 / 2 , and the difference comes about precisely because 
of the permutation signs in the expansion of dettA- a .B-fj 
[compare Eqs. (|H) and Q)]. 

Finally, the defect pair energy for the optimized posi- 
tions is simply t op t.dcf./T d = -($ max -^mm), and we see 
that Z N -i/Z N ~ £~ 1/2 exp(-e o p t .dcf./7d); the factor 
L~V 2 here represents the entropic cost for introducing 
two essentially fixed monomers into the dimer system. 

Putting everything together, we may thus write E m - ln rs 
(num. fact.) x Lexp(C op t.dof./2d). Thus, the dimer RG 
expression (E3) underestimates the smallest positive en- 
ergy by a factor L" 1 coming from the fermionic sign and 
dimcr entropic origins. 

Returning to the general case, we expect the interme- 
diate bounds, as well as the dimcr RG estimate from the 
optimized defects, to underestimate the smallest positive 
energy by some similar L~ r factors. In our numerical 
tests (see below) , we indeed find that the dimer RG pro- 
vides a strong lower bound for -Emm- However, we can- 
not be more precise in the general case because of the 
following difficulty: Unlike the random vector potential 
case, we cannot disentangle the defect energetics from 
the dimer configurational combinatorics. More explicitly, 
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the results of Rcf. [LJ suggest that a fixed pair of defects 
will have an average energy 0(+ In \r a — rp\) because the 
system constrained by the presence of the defects cannot 
take as much advantage of the low energy bonds, and 
Copt.def. will have a similar 0(+lnL) contribution in ad- 
dition to the dominant 0( — (ln_L) 3 ' 2 ) term coming from 



V. NUMERICAL TESTS 

We have performed extensive numerical tests using ex- 
act diagonalization and numerical RG. The results are 
consistent with the generalized-Gade scaling Eq. (g) in 
the delocalized phase. Although direct diagonalization 
methods do not allow us to unambiguously distinguish 
the different exponents x in the generalized-Gade scal- 
ing, indirect methods provide strong indications in favor 
of x = 3/2, different from the original Gade prediction 
of X = 2. In the localized phase, the numerical results 
are more conclusive, and all methods clearly point to a 
power-law density of states (with nonuniversal dynamical 
exponent) generated by the "string" Griffiths mechanism 
discussed earlier. 



A. Illustration 

We first demonstrate some difficulties that plague di- 
rect diagonalization studies. In this context, it is impor- 
tant to note that a strong Gade-like divergence in the 
low-energy density of states-ias never been observed in 
previous numerical studies ES'EZI In our direct diagonal- 
ization studies, we typically see power-law diverging or 
even vanishing density of states with nonuniversal expo- 
nents. This Js illustrated in Fig. where we show the 
calculatedcao density of states p(E) for the honeycomb 
lattice with hopping amplitudes chosen uniformly from 
[1 - W, 1 + W]. For W < 1, on the length scales studied, 
the density of states is apparently power-law vanishing, 
with the observed effective exponent z e g that depends 
on the strength of disorder W. It is simply impossible 
from the data to guess that there may be something else 
happening at still lower energy scales! For W = 1.0 the 
density of states is almost constant, with a slight hint on 
some weak divergence at still lower energies observed for 
the largest system studied. For W > 1.0 (the honeycomb 
lattice is still statistically isotropic in spite of the nega- 
tive hopping amplitudes allowed, but the direct dimcr 
connection is lost in this case), a power-law divergence 
is observed, and in this case, we also see some curvature 
in the log-log plot of p(E) vs E suggestive of an even 
stronger divergence at still lower energy; however, from 
these numerical results it is not possible to say anything 
more about the nature of the true low-energy singularity 
in the DOS. 




FIG. 4. Illustration of the "apparently nonuniversal" be- 
havior of the density states (normalized per lattice site) . The 
system is a 192 x 192 honeycomb lattice with hopping elements 
chosen randomly from a uniform distribution in [1 — W, 1 + W\; 
free boundary conditions are used in one direction and peri- 
odic be in the other [see Fig. |g(a)]. The data is averaged over 
10 disorder realizations; the lowest energy shown in each case 
corresponds roughly to the second-lowest energy state in the 
finite sample. Similar data for a smaller 96 x 96 system is also 
shown; from such comparisons, we conclude that the density 
of states plotted is indeed the bulk density of states. The 
effective dynamical exponent z e g is obtained by fitting the 
integrated density of states Ne to the form Ne ~ E d ' z over 
the range 2 < L 2 Ne < 200 in each case. Direction-averaged 
Lyapunov spectrum density IPLyap.o is calculated from numer- 
ical transfer matrix studies and provides a rough idea on the 
conducting properties of the system. 



A partial resolution of this discrepancy between direct 
numerics and theory lies in the fact that the asymptotic 
regime, in which the form Eq. (Eh is expected to hold, 
happens to be quite inaccessible for the exact diagonal- 
ization studies, given the computational restriction on 
the system size. Taking the original results of Ref. [I] at 
face value for the moment, we see that .the asymptotic 
form is valid only below a crossover scalecJ 



E c 



iloe 



(45) 



where Qq is some bare energy, a is the dimensionless 
conductivity of the 2D system, and y is some dimension- 
less parameter that is greater than one. The conductiv- 
ity a that enters here can be estimated as a ~ ^LyapjO 
with some 0{\) prefactor; here O^Lyap^ is the relevant 
Lyapunov spectrum density that can be obtained from a 
numerical transfer matrix analysis (see the next subsec- 
tion). In Fig. 0, we have therefore also listed the corre- 
sponding values of O^Lyap^ obtained for our system. The 
important thing to notice is that for the typically stud- 
ied systems, the estimated a is of order one and can vary 
somewhat. To get some feeling for the numbers involved, 
we set y = 1 and check how E cross changes as we vary 
W. For the three values a = 1.0, 0.5, and 0.2, we get 
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E cmss /fl = 1CT 69 , 1CT 17 , and 1CT 3 respectively! For 
energies above the crossover energy, Gade predicts an 
effective power-law density of states, and this is proba- 
bly all that has been observed in the previous numerical 
studies (but see the very recent Ref. EHl). 

Furthermore, the modified-Gade form Eq. (g) with 
x = 3/2 proposed here is a weaker divergence than the 
earlier prediction of Gade. This is probably another rea- 
son why it has proved difficult to observe anything other 
than power-law singularities in the direct diagonalization 
studies. Of course, we expect the afitual crossovers to be 
correspondingly more complicated^ but estimates sim- 
ilar to the ones obtained from Gade's original work are 
still expected to hold. The upshot of all this is that one 
should therefore be very cautious when looking for the 
asymptotic behavior in two dimensions. 

Finally, note that the foregoing suggests a way to bring 
the crossover energy up into the range accessible to direct 
numerical studies — the idea is to make the system look 
less delocalized, and this is what we turn to shortly. 



B. The choice of the system 

In the following, we present our numerical results only 
for one specific system — the brickwall (honeycomb) lat- 
tice of Fig. ||(a). Nearest-neighbor hopping amplitudes t e 
are taken randomly and independently from uniform dis- 
tribution [0, J e ). When J e = 1 for all bonds of the lattice 
(i.e., all i e 's are taken from the same distribution and the 
underlying honeycomb lattice is statistically isotropic) 
the system is not exponentially localized at E = 0, as 
can be checked by direct transfer matrix analysis. By 
allowing the horizontal bonds to be stronger on average, 
J e = e s for the dark thick bonds in Fig. ||(a) (corre- 
sponding to an anisotropic honeycomb lattice), we can 
drive the system into a localized state for strong enough 

5 > 6 C , but the system remains "delocalized" for weak 

6 < S c [see Fig. ||(b) for the phase diagram]. 

The rationale behind our specific choice of system is 
as follows: We want to perform all numerics on the same 
system, so that we can compare results of different ap- 
proaches. It is therefore best to use a system that has 
the precise connection dcttAB = Zd[A, B] with the cor- 
responding random dimer problem, so that we can fur- 
ther check our arguments of Sec. IV B. Moreover, we 



a) 



xe' 



want a system whose statistically isotropic version (i.e., 
without any enforced dimerization in values of the hop- 
ping amplitudes) shows an appreciable singularity in the 
low-energy DOS for accessible system sizes [this last re- 
quirement rules out the 7r-flux lattice — the largest z e a 
that we can observe in this system with order one bare 
randomness is z e s(L=192) = 1.5]. 

Of course, the fact that our lattice has the precise 
dimer connection does not mean that the results pre- 
sented below are non-generic. We also studied the hon- 
eycomb and square lattices with random bonds of any 
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FIG. 5. a) The brickwall lattice system we study numeri- 
cally. Different groups of hopping amplitudes are chosen from 
different uniform distributions [0, J], with J = 1 for thin lines 



(regular bonds) and J 



for dark thick lines (dimerized 



bonds), b) Phase diagram of this system from transfer ma- 
trix analysis; 8 C ~ 1.432. 



sign and obtained qualitatively similar results. In these 
systems with disorder strength of order one, we actually 
observe a stronger effective divergence in p(E), which can 
be more readily checked against the precise functional 
form, but only through multi-parameter fits, which are 
not entirely convincing. Since we can alternatively bring 
the crossover energy up into the accessible range by in- 
creasing S and simultaneously preserving the dimer con- 
nection, this is what we choose to do, in part because 
the more precise "dimer" arguments provide us with the 
strongest evidence in favor of the proposed modified- 
Gade scaling. However, we do expect that the suggestive 
strong-randomness association of the low-energy states 
in the hopping problem with the optimized defects in the 
dimer problem is likely still valid in some effective sense 
even when the fcrmionic problem does not have the pre- 
cise E = dimer connection. 



C. Transfer matrix analysis 

We now summarize the transfer matrix analysis that 
gives us the phase diagram Fig. ||(b), and also explain our 
procedure for estimating . the . conductivity a. The gen- 
eral setting is as follows :t2lcia Consider transferring the 
wave function along a strip of transverse size N = L±. 
The transfer matrix is a 2N x 2N matrix, and we calcu- 
late the Lyapunov spectrum of the corresponding random 
matrix product. The 2N Lyapunov exponents z/j come 
in ± pairs, and we focus on the positive half of the Lya- 
punov spectrum v\ > . . . > v^ > 0. The largest localiza- 
tion length £ max in the strip is obtained as the inverse 
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of the smallest (in absolute value) Lyapunov exponent, 
£ max = \jvjq. If £ max (A0 is found to increase linearly 
with N, the system is said to be delocalized (more pre- 
cisely, the 2D system is not exponentially localized). If, 
on the other hand, £ m ax(N — ► oo) saturates to a_finitc 
value, the 2D system is localized. Alternatively!^ the 
system is delocalized if there is a finite density of Lya- 
punov exponents 



? Lyap H = lim A^ 1 V 5{v - Ui) 

N^oq *• — 4 



(46) 



at v = 0; similarly the system is localized if there is a 
gap in the Lyapunov spectrum around v = 0. Moreover, 
in the delocalized state, the Lyapunov spectrum density 
^Lyap.o^ 5'Lyap(' y =0) at v = gives us a measure of the 
conducting properties of the system. 

In our bipartite hopping problem at E = 0, we ob- 
serve that the two sublattices decouple, and it is ad- 
vantageous to transfer on one sublattice only (say, the 
A-sublattice), recovering the Lyapunov spectrum of the 
whole system from appropriate symmetries. Further- 
more, for the brickwall lattice of Fig. 0(a), the sublat- 
tice transfer in the x direction involves only one column 
of A-sites at a time, but now the individual sublattice 
Lyapunov spectrum is not symmetric around zero. (In 
passing we note that in such a system with a free left 
boundary, there will be edge states living on the A-sites 
of the boundary, and the number of such edge states is 
precisely the number of negative Lyapunov exponents for 
the transfer in this case.) We employ all these simplifica- 
tions in our numerical calculations; however, for the final 
results, we always quote the TLyap,o for the whole system 
including both sublattices. 

For the brickwall lattice of Fig. g(a), there are two in- 
equivalent directions x and y. However, we have checked 
that the system is either delocalized in both directions 
or is localized in both directions. We therefore expect 
the system to exhibit true 2D behavior in the delocal- 
ized phase as long as one is not too close to the local- 
ization transition. Our numerical results are given in 
the corresponding Figs. |I] and g. To estimate ^Lyap^, 
we typically use several Lyapunov exponents that are 
closest in absolute value to zero, and we quote only the 



direction-averaged ? Lyap ,o = \J ^Lyap.o ^Lyap^o • ( If the 
Lyapunov exponents are defined as growth exponents per 
unit length instead of per lattice unit, and if the corre- 
sponding density is also defined per transverse length, the 
lattice spacings a x and a y enter both JLyap.o ~ o-x/ay 

and ^Lyap.O ~ a v/ a xi but not ?Lyap,0-) 

We find that the system is indeed delocalized for 5 = 0. 
The fact that it remains delocalized for some range of 
values 6 > can also be seen from the data at 5 = 
by considering the transfer on the ^-sublattice in the x 
direction: In this case, multiplying the horizontal bonds 
by e simply shifts the whole Lyapunov spectrum rigidly 
down by 5, Vi — > v% — 5. Thus, there remains some finite 
Lyapunov spectrum density J^yap^ for a range of 5, until 



the top of the spectrum reaches v = 0, and only then does 
the system localize. The critical 5 C can be estimated 
very accurately, since we need to know only the largest 
Lyapunov exponent in the 5 = case; for the particular 
system studied, we find <5 C « 1.432. 

What "conducting" property of the system does the. 
finite TLyap.o correspond to? Chalker and Bernhardt^ 
suggest that T^ap,!) gives, up to some numerical factor, 
the conductivity of the 2D system in units of e 2 /h : 



a 



y 



Lyap,0- 



(47) 



However, we have to ask ourselves which conductivity 
is actually being "measured" by Eq. (J47]), particularly 
since we are dealing with systems that have a singu- 
lar DOS at the Fermi level. In this respect, note that 
the above transfer matrix approach can be equivalently 
formulated as a recursion procedure for calculating the 
Green's function (E — HX^ 1 between the first and the last 
slices of the long stripj2j while the largest localization 
length from such transfer matrix analysis corresponds to 
the average modulus-squared of such a "propagator" be- 
tween the two slices. This transfer matrix estimate re- 
sembles the definition of conductivity used in the sigma- 
model literature (see, e.g., Refs. ^| and [34|), which is the 
conductivity we need in Eq. (|45|). To make our esti- 
mates more consistent, we fix the proportionality factor 
in Eq. (E7j) by using the exact results for some pure- 
systems. Thus, an elementary calculation for the pure 
7r-flux model gives JLyap.o — 2/7T, which should be com- 
pared with o 7 = 2 x 1/n for the two Dirac points of the 
model's continuum limit (see Ref . |4|) . Frora-this and simi- 
lar considerations for the brickwall lattice,E3 we conclude 
that the accurate relationship is simply er = TLyap.o- 

Finally, a couple of asides. The first concerns the ran- 
dom vector potential model, for which Ref. ^ predicts 
a = 2/ir for any strength of the randomness. This re- 
sult can be also obtained within the above transfer ma- 
trix approach, if we assume that we can indeed fix the 
unknown numerical factor in Eq. (Wq) to the appropri- 
ate constant value. When we add the random surface 
$(r) on top of, say, the pure 7r-flux model [see Eq. (jjl])], 
the E = transfered wave functions need to be simply 



multiplied by the appropriate e 



-*(a) 



or e 



*W. But this 



cannot change the Lyapunov spectrum of the transfer 
along the quasi- ID strip for any finite transverse size of 
the strip, because <I>(r) in such geometry is effectively a 
one-dimensional Gaussian surface and cannot fluctuate 
stronger than [$(iii) — $(0)] ~ Lii along the strip. At 
this stage, it is not clear to us whether a should have 
any signature of the "freezing" transition at g c observed 
in the static and dynamic properties of this system; the 
above argument seems to suggest that there is none, but 
we cannot rule out the possibility that taking the limit 
Lii — ► oo while keeping L± fixed "loses" the information 
about the 2D system. 

The second remark concerns the brickwall lattice of 
Fig. ||(a) and the transfer in the x direction. In this case, 
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FIG. 6. Density states for the anisotropic brickwall lattice 
of Fig. fcrt calculated for 192 x 192 system with open be in the 
y direction and periodic be in the x direction. The regular 
(vertical) hopping amplitudes are chosen randomly from a 
uniform distribution in [0, 1], while the dimerized (horizontal) 
hopping amplitudes are chosen from [0, e s ]. The more detailed 
characterization of the data is similar to that of Fig. 0. 



there is a direct correspondence between the conducting 
properties of the system as measured by the Lyapunov 
spectrum and the energetics of the domain walls in the 
corresponding dimer problem. Such domain walls are 
defined relative to the fully locked state that obtains for 
5 > 8 C , and are forced to run in the x direction and can- 
not terminate; their number at a given 6 is precisely the 
number of positive Lyapunov exponents in the hopping 
problem, while the free energy cost per unit length of 
adding or removing domain walls is precisely the spacing 
between the Lyapunov exponents near v = 0. 



D. Exact diagonalization studies: density of states 

We also monitor the density of states as we increase 
dimerization from S — 0, eventually driving the system 
into a localized state for S > S c . The corresponding re- 
sults are shown in Fig. pi together with information about 
the system's localization properties. Consistent with our 
earlier discussion, the DOS divergence is seen to become 
stronger as we approach the localization transition, with 
the effective exponent z c g(L=l92) peaking at the tran- 
sition. Deep in the localized phase, the log-log plot of 
p(E) vs E is indeed a straight line, and the quoted val- 
ues of z give the actual bulk dynamical exponents. In 
the delocalized phase, on the other hand, there is still 
some curvature at the lowest energy scale, and the ex- 
tracted values z e ff can only serve as rough indicators of 
the strength of the divergence. 

The actual data analysis is performed as follows: In the 
delocalized phase, we fit the integrated density of states 
Ne to the generalized-Gade form 



N E = aT' e ® 



(48) 



where Te = \o{£Iq/E). Specifically, we try two values for 
the exponent x — the original Gade x — 2 and the mod- 
ified x = 3/2. The included prefactor Tj x is based on 
an assumption that there are no additional corrections 
to the asymptotic E — > form p[E) ~ E~* exp(— Tj x ). 
This is something that we do not really know — in fact, 
our lowest-gap analysis of Sec. IV B| seems to suggest that 
there may be corrections stronger than any such T r pref- 
actor. In this situation, we treat the above fit function 
only as a baseline. To be consistent, we should also allow 
for the uncertainty Te — ► r^+Fo in the bare energy scale 
relative to which the log-scale T is defined. As one might 
suspect, for a fixed fit range, we can approximate any 
DOS curve with the three parameters a, c, and To, so this 
sort of analysis is not really conclusive. What we can do 
however is to look for the overall consistency of such fits: 
We can ask, for example, how the fit parameters change 
as we change the fit region, or how sensible the obtained 
parameters are (thus, the parameter a should be of order 
one). Generally, the x = 3/2 form fares somewhat better 
in such analysis; also, the extracted numerical value of 
the parameter c from such fits compares favorably with 
the lowest-gap results of the next subsection. But this 
is as far as our direct DOS studies can take us in the 
delocalized phase. 

In the localized phase, we have also tried fitting the 
DOS using a power-law times a logarithmic correction 
N E ~ E d l z (T E ) r - We typically find that the best fit cor- 
responds to r = 0, suggesting that the density of states 
has a simple power-law form in the localized phase. This 
is also our conclusion from the lowest-gap studies pre- 
sented below. 



E. The lowest energy state: exact diagonalization 
and dimer optimized defects studies 

We can extend our numerical tests further by consid- 
ering the distribution of the smallest positive energy in 
finite samples, paying particular attention to scaling with 
sample size. Such a direct diagonalization study can be 
performed for a factor of two larger systems than in our 
DOS studies, while indirect dimer methods can take us 
another factor of four in system size. This well-controlled 
and rather sensitive numerical approach als o corr esponds 
closely to our analytical arguments in Sec. IV B| . 

We consider the same brickwall lattice system of 
Fig. ra(a); however, we now consider (2M X + 1) X 2M y 
"odd x even" samples with spiral boundary conditions 
in the x direction and free boundaries in the y direc- 
tion. The specific choice of the boundary conditions 
is such that the dimer equivalence det £ab = Zd [A, B] 
holds precisely for this system;E3 this allows us to com- 
pare the lowest-energy state of the hopping problem with 
the optimized defects in the corresponding di mer p rob- 

discussion of Sec. IV B . To 
3 the dimer RG "estimate" 



lem, exactly paralleling tt 
this end, we also calculated 
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7. "Lowest-gap" analysis in the delocalized phase 
of the system of Fig. ]A We consider (L + l)xL finite 
samples with spiral be in the x direction and free be in the 
y direction. Open symbols show distributions of the smallest 
positive energy calculated by exact diagonalization methods, 
while filled symbols show distributions of the corresponding 
dimer RG "estimate" Eq. 
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FIG. 8. Detailed analysis of the distributions of Fig. [71 The 
mean r opt (L) and the standard deviation cr(r op t; L) are plot- 
ted with open and filled symbols for the exact diagonaliza- 
tion and dimer RG results correspondingly. Note the loga- 
rithmic scale for both the In L and r opt axes. To gauge the 
L-dependence of F op t and a(r op t), we also plot different pow- 
ers (InL)* with x = 2, 3/2, 1, and 1/2. 



Eq. (E3I) for the lo west g ap, which is the quantity used in 
the bounds of Sec. IVB. 



For the purposes of illustration, we show our results 
for two values of 8, one in the delocalized phase (6 = 0.8) 
and one in the localized phase (8 = 2.0), both chosen well 
away from the transition point (<5 C = 1.432). 



1. Delocalized phase, 5 = 0.8 

Fig. shows distributions of the logarithm of the 
smallest positive energy in our samples in the delocal- 
ized phase. The exact diagonalization results for r op t = 
ln(l/L7+ in ) are shown with open symbols and extend up 
to system sizes L = 384, while the corresponding dimer 
RG estimates — C op t.dcf. [cf. Eq. (p3|)] are shown with 
filled symbols and extend up to L = 1536. We see that 
the distributions move out strongly to larger absolute 
log-energies (i.e., smaller energies) with increasing L and 
broaden somewhat at the same time. We also see that 
the dimer RG estimates for the lowest gap are consis- 
tently smaller and indeed provide strong lower bounds, 
at least for the larger system s. This is in agreement with 
our discussion in Sec. IVB. 

In Fig. @, we analyze the behavior of the distributions 
P(F opt ;L) of Fig. 0, plotting their mean r opt (L) and 
the standard deviation cr(r opt ; L). On a plot with linear 
scales for both InL and r opt (not shown), we observe a 
visible curvature r opt ~ (\x\L) x with x > 1, indicative of 
a stronger than any finite- z singularity In Fig. |8j, we use 
logarithmic scales for both In L and r opt to get a rough 
estimate of the exponent x. The exact diagonalization 
results are seen to fall between the x = 1 and x = 3/2 
lines, well away from the Gade x = 2; at the same time, 



the optimized defect results, which extend to larger sys- 
tem sizes, clearly approach the x = 3/2 line, as has been 
already checked in Ref. [lj for a related model. We also 
note a much weaker ~ (InL) 1 / 2 dependence of the width 
of the distributions on the system size, consistent with 
the results of Ref. Il2| (summarized in Sec. IVB). 

In a more quantitative analysis, following Ref. |l 
perform linear fits for [r opt ] 2 / 3 and [<r(r opt )] 
tions of In L . For the exact diagonalization data, we 
obtain the asymptotic scaling r opt « 0.81 (InL) 3 / 2 and 
a sw 0.42 (InL) 1 / 2 . We can compare this with the result 
from the DOS studies: performing the generalized-Gade 
fit Eq. (Q) of the data of Fig. ra, we extract the param- 
eter c = 2.2, which translates to a fairly close asymp- 



we 
2 as func- 



totic scaling prediction r opt 



0.87 (InL) 3 / 2 . Similarly, 
the dimer RG data yields r opt ~ 1.30 (InL) 3 / 2 , not un- 
reasonably far from the above exact diagonalization pre- 
dictions [as discussed in Sec. IVB, we expect stronger 



subleading O(lnL) terms for the exact lowest gap r opt 
than for the dimer — C op t.dcf.]- Finally, we note that 
Ref. [39] calculated bulk defect density in a related vor- 
tex glass model as a function of the defect core energy 
(which corresponds to r opt in the fermion problem), ob- 
taining Nr ~ exp (— cr 074 ), fairly close to the expected 
Nr ~ exp (—cT 1 ^ x ) scaling form with x = 3/2. 



2. Localized phase, S = 2.0 

The lowest-gap analysis is simpler in the localized 
phase, and is shown in Figs. |9-10, In this case too, the 
distributions move out to larger r opt with increasing L, 
but now do so at a slower pace. We clearly observe sim- 



15 





5=2.0 
Lowest Gap 


1 












diag 




dimer 




0.4 


L=24 -e 
96 -e- 


■ 

• 


24 
96 






384 ^ 




A 
▼ 


384 
1536 












03 


■, ^ 






. 




iA p-n *.;\ t--t 






fm ;<?' \<s - /*\ ; '■• 






1 i // Us, i: \\ ; 




0.2 


1 / *\ ! i *\ ; t 
r vv •■••■ ?■'. 




0.1 
00 


// Ma i& *'•■. \ 





10 15 

r 0Dt = ln(l/Et 



20 



25 



- opt ~ '" v " '-mm' 
FIG. 9. Lowest-gap analysis in the localized phase, 8 = 2.0, 
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pie finite z scaling r opt = const + z InL , with z = 2.56, 
which compares well with the estimate of z — 2.64 from 
the DOS studies. Also, the widths of the distributions 
essentially do not change with system size. Note that 
the effects discussed in Sec. IV B produce a much smaller- 
difference between the exact diagonalization results and 
the dimer RG estimates in the localized phase, since the 
distance between the optimal defects is expected to scale 
very weakly ~ In L with the sample size. Thus, at least 
at S — 2.0 deep in the localized phase, and for the sizes 
studied, the dimer RG gives numerically accurate results 
for the lowest gap and the two peaks of the corresponding 
wave function in each sample. 



VI. DISCUSSION 

Our main results have already been summarized in the 
Introduction, and here, we confine ourselves to some un- 
resolved questions in these bipartite models, as well as 
some issues in the closely related ImRH problems. 



Our first remark concerns the zero-energy wave func- 
tion in the "delocalized" critical phase of the general bi- 
partite problem. The multifractal properties of this wave 
function are expected to be trivial r(q) = — this is be- 
cause of the flow to strong disorder g — > oo (see, e.g., 
Ref. [L0|). Equivalently, since the logarithm of the wave 
function magnitude is essentially the potential v(r) seen 
by the defects of the corresponding dimer problem [see 
Eq. (39)], none of the participation ratios depends on 
system size — this is a consequence of the ground state 
dominance for this random potential. In this respect, 
the critical wave function looks fairly localized. The typ- 
ical correlations like Ct yp (r) = exp 



log|*(r)tf(0)| 



are 

thus directly related to the correlations of v(r). On the 
other hand, the situation is less clear when we consider 
average correlations like C av (r) = |$(r)$(0)|. Estimat- 
ing these requires a better understanding of the statistics 
of such random surfaces. For instance, it is not clear if 
Cav (f) is dominated by the near- returns of the surface to 
its global extrema, as is the case in the one-dimensional 
chain. Moreover, if this is indeed the case, what is the 
probability of finding such quasi-degenerate minima at a 
fixed distance r? 

Another feature of the delocalized phase that we have 
not discussed is the suggestedBa universal susceptibility 
variations in the vortex glass problem, and it would be 
interesting to study the corresponding response in the 
fcrmionic problem in more detail. 

The character of the transition between the delocalized 
and localized phases is a completely separate issue that 
we have ignored altogether. It is not obvious if this crit- 
ical end-point of the line of fixed points that character- 
ize the metallic phase exhibits dynamical scaling distinct 
from the metallic phase. Another interesting feature of 
this transition is that the conductivity a xx remains fi- 
nite (because the Lyapunov spectrum density is finite at 
the top of the Lyapunov band), while ~5 yy vanishes. It is 
not clear how this extreme anisotropy affects the proper- 
ties of the transition, and whether anything interesting 
remains. 

In the localized phase, we have not discussed in de- 
tail the structure of our low-energy "string" excitations. 
Presumably, they look similar to the domain walls (flux 
lines) in the corresponding random dimer (vortex glass) 
problem. Note that in the particular localized phase that 
we studied, the strings are forced to run in the x direc- 
tion. In a more general localized phase, obtained by in- 
troducing some other dimerization pattern, the situation 
is more complicated, but it still seems that it will be 
some kind of directed strings that will contribute most 
to the low-energy density of states: The strings need to 
be stretched so that the end-to-end distance (setting the 
tunneling frequency scale) is of order the string length 
(determining the occu rrenc e probability) for our count- 
ing arguments of Sec. [II B| to work. An interesting sys- 
tem to study in this respect is a 2D spinless supercon- 
ductor with time-reversal invariance, which maps onto 
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a bilayer BPRH. When the "p^-wave" superconducting 
pairing amplitude Ajj+x is zero, the system is just a 
doubled Anderson localization problem, whereas increas- 
ing A from zero, we obtain a bipartite localized phase 
"connected" to the standard Anderson insulator point. 
Note the richness of the full phase diagram in this sys- 
tem: a simple Lyapunov spectrum shift argument shows 
that increasing A further drives the system from this lo- 
calized phase first into a critical delocalized phase, and 
then again into a localized phase similar to our staggered 
band insulator! 

We conclude with some comments regarding more 
generic ImRH problems. As is clear from the preced- 
ing discussion, we expect similar Griffiths string effects 
in the localized phases of such systems. An important 
question is, of course, whether this Griffiths mechanism 
can compete with other mechanisms of "filling the gap" 
(which can, for example, produce a constant contribu- 
tion to the density of states at the band center as in 
Rcf. ph. In particular, are there situations with a power- 
law divergent density of states that has its origins in such 
Griffiths effects? While we have not studied these ques- 
tions in detail, there are certainly situations where this 
does indeed happen. .-Thus, in the particular fermionic 
ImRH representation^ of the two-dimensional random 
bond Ising model, our preliminary results suggest that 
the density of states is power-law vanishing to the left 
(the less disordered side) of the Nishimori line and power- 
law diverging to the right (this particular example was 
suggested to us by Nick Read and the corresponding re- 
sult proved for a one-dimensional toy- model in Ref. 42). 



Finally, 2D ImRH systems are believed to also have true 
metallic phases,Q and it would be interesting to consider 
these from this perspective. In particular, can one use an 
analogous dimer connection to gain further insight into 
their properties? 
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APPENDIX A: 2D GAUSSIAN SURFACE 

The following facts are used in the main text to ob- 
tain bounds for the dynamical scaling in the random vec- 
tor potential model; these are transcribed directly from 
Rcf. 11 for the two-dimensional Gaussian surface $(r), 
Eq. (15), in an L x L box (stipulating fdr$(r) = 0). 

Extremal properties of the surface <&(]") are sharply 
defined: e.g., the maximum scales as 



$ max (L) = 2 



InL 




with the sample-to-sample dependence entering only 
through an O(l) x ^fg random variable t><I> m ax- 
The following "partition function" 



3 = E« 



-2#(r) 



has a sharply defined logarithm: 
ln3 9 (L)=2(l + -^-)lnL + A g , 

ln3 9 (£) =41n£- -lnlnZ + A, 



(A2) 

g < g c = 2ir ; 
9 = 9c ; 



ln3 9 (£) = J— (4rnL--lnlnZ,) + A ff , g > g c ; 
V 9c & 

where A ff is an O(l) random number (different in each 
case). Note that for strong disorder g > g c , J>{L) is dom- 
inated by the global minimum: 



3,(1) = e 



_ „-2* mill (L) 



£■ 



-2[*(r)-* min (L)] 



(A3) 



where the sum contributes no L-dcpcndcnce because 
e -2[$(r)-$ min ] remams normalizable in the limit L — > oo. 



In In L + 5$„ 
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